%% Weighing of 5 objets
%%
% A cook wants to evaluate the effects of 5 factors of the recipe of the Opera cake. Each factor can be varied in a specific range without risking to alter the process and lost a part of the production. The effects are measured by a panel of testers. The list of factors
% Quantity of cream
% Quantity of ganache of chocolate
% Cooking time for the biscuit
% Quantity of sirup for punching
% Alcohol
% Different strategies are evaluated with the help of the DOE theory.
% Table of Contents
% Weighing of 5 objets
%
%     Strategy 1 (offset, and 1x1)
%     Strategy 2 (all, then 1x1)
%     Strategy 3 (Offset, then 2x2)
%     Strategy 4 (All, then 2x2)
%     Strategy 5 (Offset, then 3x3)
%     Strategy 6 (All, then 3x3)
%     Strategy 7 ( Hadamard)
%     Evolution of the trace
%     Evolution of the determinant
%     Conclusion
%     Function of design analysis
%
%
%% Strategy 1 (offset, and 1x1)
% This strategy is called sometimes OFAT and is known to be in efficient.
%
% Matrix of experiments E

E1=[0 0 0 0 0;eye(5)];
disp(E1)

%
% Design evaluation

% call the function 'analyse' that normalize the design, compute the
% trace, the determinant and the dispersion matrix
[T(1),DET(1),Dii(1),VIF(1)]=analyse(E1,'Strategy 1: Offset, then 1x1');

% Observation on strategy 1
% We observe that the variance on the offset a_o is multiplied by 3.5 in comparison with the experimental variance.
% The viariance of the other coefficients of the model is equal to half of the experimental variance, whereas the optimun is 1/6.
%
%% Strategy 2 (all, then 1x1)
%
% Matrix of experiments

E2=[1 1 1 1 1;eye(5)];
disp(E2)

%
% Design evaluation

% call the function 'analyse' that normalize the design, compute the
% trace, the determinant and the dispersion matrix
[T(2),DET(2),Dii(2),VIF(2)]=analyse(E2,'Strategy 2 : All, and 1x1');

% Observations on strategy 2
% The variance of the model coefficients corresponds to 22% of the experimental variance, while the optimun is 17% (1/6).
%
%
%% Strategy 3 (Offset, then 2x2)
%
% Matrix of experiments

E3=[0 0 0 0 0
   1 1 0 0 0
   0 1 1 0 0
   0 0 1 1 0
   0 0 0 1 1
   1 0 0 0 1];
disp(E3)

%
% Design evaluation

% call the function 'analyse' that normalize the design, compute the
% trace, the determinant and the dispersion matrix

[T(3),DET(3),Dii(3),VIF(3)]=analyse(E3,'Strategy 3: Offset, then 2x2');

% Observations on strategy 3
% The variance of the model coefficients corresponds to 22% of the experimental variance, while the optimun is 17% (1/6).
%
%% Strategy 4 (All, then 2x2)
%
% Matrix of experiments

E4=[1 1 1 1 1
   1 1 0 0 0
   0 1 1 0 0
   0 0 1 1 0
   0 0 0 1 1
   1 0 0 0 1];
disp(E4)

%
% Design evaluation


[T(4),DET(4),Dii(4),VIF(4)]=analyse(E4,'Strategy 4: all, then 2x2');

% Observations on strategy 4
% The variance of the constant corresponds to 17% of the experimental variance, this is the optimum.
% The variance of the other coefficients corresponds to 33% of the experimental variance, which is far from the optimum of 17%.
%
%% Strategy 5 (Offset, then 3x3)
%
% Matrix of experiments

E5=[0 0 0 0 0
   1 1 1 0 0
   0 1 1 1 0
   0 0 1 1 1
   1 0 0 1 1
   1 1 0 0 1];
disp(E5)

%
% Design evaluation

[T(5),DET(5),Dii(5),VIF(5)]=analyse(E5,'Strategy 5: Offset, then 3x3');

% Observations on strategy 5
% Idem strategy 4.
%
%% Strategy 6 (All, then 3x3)
%
% Matrix of experiments

E6=[1 1 1 1 1
   1 1 1 0 0
   0 1 1 1 0
   0 0 1 1 1
   1 0 0 1 1
   1 1 0 0 1];
disp(E6)

%
% Design evaluation

% call of the function analyse
[T(6),DET(6),Dii(6),VIF(6)]=analyse(E6,'Strategy 6: Offset, then 3x3');

% Analysis of strategy 6
% The variance of the model coefficients corresponds to 37.5% of the experimental variance, while the optimun is 1/6 \approx17 \%
%
%% Strategy 7 ( Hadamard)
%
% In this strategy we make use of the Hadamard matrix. As it does not exist an Hadamard matrix for 5 factors, we generate the matrix for 8 factors and use only the 6 first columns.

% create a matrix of Hadamard 8x8
E7=hadamard(8);

% select the columns 2 to 6 for coding the factors
E7=E7(:,2:6);
disp(E7)

%
% Design evaluation

% call of the function analyse
[T(7),DET(7),Dii(7),VIF(7)]=analyse(E7,'Strategy 7: Hadamard');

% Analysis of strategy 7
% Optimal design: variance of the coefficients corresponds to 1/8 of the experimental variance for 8 experiments.
%
%% Evolution of the trace

figure
plot(T,'ok','MarkerSize',12,'MarkerFaceColor','r')
title('Evolution of the trace', 'FontSize',16)
ylabel('Trace of the matrix of dispersion','FontSize',14)
xlabel('Strategies','FontSize',14)
grid on
xticklabels({'offset, 1x1','all, 1x1', 'offset, 2x2','all, 2x2','offset, 3x3', 'all, 3x3', 'Hadamard'})

%
%% Evolution of the determinant

figure
plot(DET,'dk','MarkerSize',12,'MarkerFaceColor','b')
title('Determinant evolution', 'FontSize',16)
ylabel('Dispersion determinant','FontSize',14)
xlabel('Strategies','FontSize',14)
grid on

xticklabels({'offset, 1x1','all, 1x1', 'offset, 2x2','all, 2x2','offset, 3x3', 'all, 3x3', 'Hadamard'})

%
%% Conclusion
%

S=table({'offset, 1x1';'all, 1x1'; 'offset, 2x2';...
    'all, 2x2';'offset, 3x3'; 'all, 3x3'; 'Hadamard'},...
    'VariableNames',{'Strategy'});
S.Strategy=categorical(S.Strategy);
S.Runs=[6 6 6 6 6 6 8]';
S.Trace=T';
S.Dii=Dii';
S.VIF=VIF';

% The table here below summerize the characteristics of each strategy.
% 'Trace' correspond to the trace of the dispersion matrix
% 'Dii' corresponds to the value of the diagonal element of the dispersion matrix of the effects (not the constant).
% 'VIF' corresponds to the variance inflation factor of the effects (not the constant).

disp(S)

% Strategy 2 (all and 1x1) has the best performances with 6 runs.
% Strategy 7 (Hadamard) has the best performances but 8 runs.
%
%
%% Function of design analysis

function [T,DET,Dii,VIFi]=analyse(E,title_string)

% Axis characteristics
Axe1=[0 7 0 4];

% Normalisation of the matrix of experiment
E_st=rescale(E,-1,1,"InputMin",min(E),"InputMax",max(E));

disp('Normalisation of the matrix of experiment E ')
disp(' ')
disp(E_st)

% Compute the model matrix
X=x2fx(E_st,'linear');

% Compute the matrix of dispersion
D=inv(X'*X);
disp('Matrix of dispersion D')
disp(' ')
format bank
disp(D)
format short

figure
imagesc(D)
colorbar

% barchart of the diagonal elements
figure
bar(diag(D))
title(title_string, 'FontSize',16)
ylabel('Diagonal elements of the matrix of dispersion','FontSize',14)
xlabel('Coefficients','FontSize',14)
axis(Axe1);

% VIF
VIF=diag(inv(corrcov(D)));
disp('Variance inflation factors')
disp(' ')
disp(VIF)
VIFi=VIF(6);



% Compute the trace and the déterminant of the dispersion matrix
T=trace(D);
DET=det(D);
Dii=D(5,5);
end